Kondo effect and channel mixing in oscillating molecules 
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We investigate the electronic transport through a molecule in the Kondo regime. The tunnel- 
ing between the electrode and the molecule is asymmetrically modulated by the oscillations of the 
molecule, i.e., if the molecule gets closer to one of the electrodes the tunneling to that electrode 
will increase while for the other electrode it will decrease. We describe the system by a two-channel 
Anderson model with phonon-assisted hybridization. The model is solved with the Wilson numer- 
ical renormalization group method. We present results for several functional forms of tunneling 
modulation. For a linearized modulation the Kondo screening of the molecular spin is caused by 
the even or odd conduction channel. At the critical value of the electron-phonon coupling an un- 
stable two-channel Kondo fixed point is found. For a realistic modulation the spin at the molecular 
orbital is Kondo screened by the even conduction channel even in the regime of strong coupling. 
A universal consequence of the electron-phonon coupling is the softening of the phonon mode and 
the related instability to perturbations that break the left-right symmetry. When the frequency of 
oscillations decreases below the magnitude of such perturbation, the molecule is abruptly attracted 
to one of the electrodes. In this regime, the Kondo temperature is enhanced and, simultaneously, 
the conductance through the molecule is suppressed. 

PACS numbers: 72.15.Qm,73.23.-b,73.22.-f 



I. INTRODUCTION 

The Kondo effect, a generic name for processes related 
to an increased scattering rate off impurities with in- 
ternal degrees of freedom, reveals itself in mesoscopic 
systems as increased conductance at biases and tem- 
peratures low compared to the Kondo temperature. It 
has been observed in measurements of transport through 
quantum dots^, atoms, and molecules^ii^ii^i^'^iii^. Specific 
to molecules is the coupling of electrons to molecular os- 
cillations. The molecular internal vibrational modes and 
oscillations of molecules with respect to the electrodes 
have been proposed to account for the side-peaks in the 
non- linear conductance*'^'^. In addition, the electron- 
phonon coupling can explain the anomalous dependence 
of the Kondo temperature on changing the gate voltage 
at zero bia o^'^'^°i^^ . Since the electrode-molecule junc- 
tions are candidates for devices such as molecular diodes, 
switches, and rectifiers the research in this field is increas- 
ing despite its complexity and the difficult experimen- 
tal characterizationi2iHiil. Recently, also the notion of 
quantum phase transition was introduced in the analysis 
of such systemsi^. We believe that for the interpretation 
of the experimental results a better understanding of the 
behavior of simple theoretical models is necessary. 

Here we study the influence of the electron-phonon 
coupling in the Kondo regime where a single molecular 
orbital is occupied on average by one electron. We con- 
centrate on the case where the oscillations of the molecule 
with respect to the electrodes affect the tunneling as de- 
picted schematically in Fig. [TJ That is, the tunneling 
toward the left and right electrodes is given by overlap 
integrals Vl,r{x) that arc modulated by the displacement 
X of the molecule from the mid-point between the elec- 
trodes asymmetrically: Vr increases and Vl decreases for 
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Figure 1: (Color online) Schematic plot of the model device. 



X positive, and opposite for x negative. 

Assuming the electrodes are identical, it is conve- 
nient to introduce symmetric and antisymmetric com- 
binations of states in the electrodes. With respect to 
inversion, they form even and odd conduction channels. 
The odd channel is coupled to the molecule only due to 
the asymmetric modulation of tunneling. For example, 
in the linear approximation Vl,r{x) = V{\ =F gx) - a 
prefactor V sets the intensity of the tunneling and the 
electron-phonon coupling constant g its variation due to 
the displacement^^ - the even channel is coupled to the 
molecule directly and the odd channel is coupled to the 
molecule via a term proportional to gx. 

As a consequence of coupling the molecular orbital to 
two channels the low-energy behavior is that of the two- 
channel Kondo (2CK) model^iii. The screening of the 
spin occurs in the channel with the larger coupling con- 
stant. If the couplings match, an overscreened, i.e., a 
genuine 2CK problem with a non-Fermi liquid behaviour 
results. For a linearized model such a fixed point has 
indeed been found with simulations based on numerical 
renormalization group^^. 

Of interest is also the renormalization of phonon fre- 
quencies. Quite generally, the characteristic frequency of 
the oscillations decreases with increasing electron-phonon 



coupling. In the Anderson-Holstein model the softening 
of the phonon mode is related to the increased charge 
susceptibilit y^^'^" , which occurs due to dynamical break- 
ing of the particle-hole symmetry for negative effective 
repulsion U . In the present case the softening occurs as 
well and it is related to the dynamical breaking of in- 
version symmetry?. Due to the softening, the instability 
towards perturbations breaking the symmetry emerges. 
On the mean field level^i, the instability is seen as an 
asymmetric ground state with large average x in systems 
with inversion symmetry. 

In this work we extend the existing analysis in two 
ways, (i) Motivated by the lack of the inversion sym- 
metry of typical experimental devices we include the in- 
version symmetry breaking perturbation, (ii) We check 
which features persist if the tunneling is taken to depend 
exponentially on the displacement: Vl,r{x) cx expi^gx). 
In particular, it is shown that the softening of the phonon 
mode and the corresponding instability occur universally 
but the 2CK fixed point appears as an artifact of the lin- 
earization only. It is shown also that the softening is due 
to the kinetic energy gained by the dynamical breaking of 
inversion symmetry and that it occurs also for vanishing 
repulsion, U ^ 0. 

Both directions of research have been pursued in the 
context of nano-electromechanical system o^^i^^'^** , where 
only the lowest orders in tunneling are considered and 
the Kondo correlations are thus lost. A similar approach 
has been followed in analyzing the influence of pair tun- 
neling for negative effective U in the Anderson-Holstein 
mode l^^i^^ . On the other hand, the influence of the expo- 
nential dependence of tunneling rates on x in the Kondo 
regime has been analyzed in Ref. l23- However, the dis- 
placement X is not treated as a dynamical variable there 
but only as an external control parameter. 

The paper is organized as follows. In the next section 
we describe in more detail the models under considera- 
tion. We have performed the numerical calculations us- 
ing Wilson's numerical renormalization group (NRG) and 
projection operator method of Schonhammer and Gun- 
narsson (SG) which we briefly describe in Section III. In 
Section IV we present analytical and in Section V numer- 
ical results. We conclude by critically commenting the 
obtained results and their applicability. A comparison 
between the NRG and SG results is given in Appendix A 
followed by Appendices B and C containing the deriva- 
tions of the Schrieffer- Wolff transformation and the con- 
ductance formulas. 



II. MODELS 

We model the system with the Hamiltonian 

H = ii/mol + Hl + Hjj + i?vib + H , (1) 

where iJmoi describes an isolated molecule, Hl and Hr 
the left and the right electrode, respectively, iJvib a vi- 
brational mode and H' the phonon-assisted coupling of 



the molecular orbital to the electrodes. The molecule 
consists of a single orbital with energy e, which is in ex- 
periment modulated by the gate voltage. The repulsion 
between two electrons simultaneously occupying the or- 
bital is U, 



Hmoi = e("-T + "i) + Un-^ni, 



(2) 



where the number operators n„ — d^d^r count the num- 
ber of electrons in the orbital with spin a =t7i- The 
symbols c^^^d'^' denote electron annihilation (creation) 
operators in the electrodes and molecular orbital, respec- 
tively. We are here interested in the particle-hole sym- 
metric point e = — C//2 only. The oscillator part is 



H. 



na)c 



(3) 



describing the oscillations with frequency Q. and a^ is 
the boson creation operator. The left and right elec- 
trodes are described by bands of noninteracting electrons 
Ha = ^ka ^knkaa for a = L,R, respectively, where 
nkatr = CkarjCkacr counts the electrons with spin a and 
wave vector k; tk is the dispersion of the band in the 
electrode a. The chemical potential is set to the mid- 
dle of the band (/i — 0) corresponding to the half-filled 
regime where the molecule is on average singly occupied. 
In NRG calculations a flat band with constant density of 
states p — l/(2£') and in SG calculations a tight-binding 
band with p{uj) = 1/ (iry/D'^ — cj^j are used, where D is 
the half- width of the band. 

The tunneling between the molecular orbital and the 
electrodes, which is described by 



H' = VLix)dL + VB.ix)vR, 



(4) 



occurs via the hybridization operators (assuming here the 
tunneling is ^-independent) 



/ .^kaa'^ 
ka 



(5) 



multiplied by the overlap integrals Va{x) = Va{a + a^), 
where the displacement is explicitly quantized. 

It is practical to define even and odd combinations of 
states in the electrodes, respectively 



/K \ykLa ^ ^fcJ 



(6) 



In this basis W reads 

H' ^V,{x)v, + Vo{x)vo (7) 

where 

VL{x)±Vn{x) 



VeA^) = 



%/2 



(8) 



modulate the tunneling to even and odd channels. Hy- 
bridization operators Vf,^o correspond to Eq. ([S]) for a = 
e, o, respectively. Note that 



|K(x)| > \Vo(x)\ 



(9) 



if Vl,r{x) are both positive or both negative for all x. 

In this paper we perform the calculations using several 
functional forms of Va {x) depicted in Fig. [2l 



A. Overlap integrals 

In a realistic experimental situation the tunneling be- 
tween the molecule and the tip of an electrode will be 
saturated at small distances and it will progressively de- 
crease with increasing distance of the molecule from the 
electrode. The precise functional dependence of overlap 
integrals will in general depend on details of the molecule 
and the tips of the electrodes, but the overall behavior 
should be as shown in Fig. [2lja) with dotted line. 



1. Linear modulation 

The simplest form of overlap integrals is obtained by 
the expansion to lowest order in displacement resulting 
in linear modulation (LM) 



VL.,R{x)^v[iT{gx + 0]^ 



(10) 



The tunneling matrix element, constant V for g = 0, is 
linearly modulated by displacement for g > 0. We as- 
sume the system is almost inversion symmetric. A small 
^ > is the magnitude of the symmetry breaking per- 
turbation. In the symmetrized basis the overlap integrals 
take on the following form 



Ve^V2V, Vo^V2Vigx + 0- 



(11) 



Note that Eq. (Ilip docs not satisfy the requirement 
Eq. © for gx > 1 — (, because the overlap to the left 
electrode becomes negative and its absolute value starts 
to increase with increasing x (dashed region in Fig. [2]). 
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Figure 2: (Color online) (a) Various forms of the tunneling- 
modulation. The unphysical regime of LM where the tunnel- 
ing starts to increase with increasing distance to the electrode 
is indicated by dashing, (b) The breakdown of EM. In the 
right-hand part the effective potential drops without bounds. 
This is regularized by the phonon cutoff, which corresponds 
to the hard-wall boundary. Two different cutoff regimes are 
indicated by dotted and dashed area. 



gx. Namely, the model with EM is unstable towards 
large displacements as can be understood by the following 
simple argument. For large x, the largest energies in the 
problem are V e-scplgx) and flx^. This limit corresponds 
to two sites coupled by a tunneling term with tunneling 
proportional to exjp{gx). The total energy of one electron 
on these two sites is E ^ ilx'^ /A — V ex-p{gx) . Therefore 
the oscillator in the EM moves in an effective potential 
of the form depicted in Fig. [Hb) (full), unbounded from 
below for large |x|. 



2. Exponential modulation 



Limitations of models 



A better approximation to the overlap integrals could 
be exponential decay at large distances and arguably 
more realistic model is given with the exponential mod- 
ulation (EM) of tunneling by 



VlM^) = V [exp(=Fffa;) =F C] , 
or equivalently 



(12) 



Ve = V2V cosh{gx), Vo = V2V [sinh(52;) -I- C] , (13) 

which are positive and for C — manifestly satisfy the 
relation Eq. ©. By expanding the couplings to lowest 
order in x, the LM is recovered. 

The EM eliminates the negative overlap but introduces 
another problem due to divergence of exp(gx) at large 



In real systems, the overlap integrals will be neither 
negative nor divergent. The large x behaviour of EM can 
be corrected by adding higher terms in displacement to 
the oscillator potential, which corresponds to hardening 
of the 'spring' for large x. In our numerical calculations 
such a hardening is incorporated in the form of a phonon 
cutoff which acts as a hard-wall boundary, thereby elim- 
inating the states corresponding to displacements larger 
than ~ 2 vZ, where L is the maximal number of phonons 
allowed. In Fig.[21^b), the dotted and the dashed regions 
indicate two different cutoff regimes. For the larger cutoff 
also the resulting effective potential is sketched (dotted). 

By incorporating the cutoff into EM the results qual- 
itatively depend on additional parameter L because the 
choice of cutoff determines the form of the effective po- 
tential near low energies. However, without some kind 



of a regularization of the model the model with EM is ill 
defined; we show later that the average displacement (or 
its fluctuations for C = 0) diverge for all 5 > at L -^ oo. 
We note that there may be several cases, where the 
modulation is not a simple function. To describe the ex- 
perimental situation, it might even be needed to include 
the anharmonicity of the potential as well. However, a 
convenient starting point is to first clarify specific regimes 
of simplified models and the consequences of the approx- 
imations. In this paper, we first analyze the model with 
LM comprehensively. Later we discuss EM for a specific 
phonon cutoff to highlight which of the results obtained 
using LM are artifacts of the linearization. Finally, the 
exponential divergence of overlap integrals is regularized. 
Fig. [2)Ja) (dotted), and it is shown which results persist 
also for this model. 



eters of the auxiliary Hamiltonian (for example, by fix- 
ing Vl/Vr = r), the minimization procedure gives states 
l^r), for which the (\I/r.|x|^,) ~ Xj. is in general finite, 
and energies E^ = {^r\H\^r)- The effective potential is 
estimated by pairs {xr,Er). 



IV. ANALYTICAL RESULTS 

We studied the model Eq. IH) numerically and the re- 
sults are presented in Section V. Nevertheless, from an- 
alytical results in special limits we anticipate different 
regimes of behavior and the values of parameters where 
these regimes emerge. 



Linear modulation 



III. NUMERICAL METHODS 

Most of the numerical results presented here have 
been obtained using the Wilson numerical renormaliza- 
tion group2^i2S(NRG) method. The NRG procedure is 
based on adding sites to the system iteratively with hop- 
ping matrix element to the nth added site decreasing as 
j^-n/2^ At each step the resulting Hamiltonian is diag- 
onalized and lowest K eigenstates are kept. The expo- 
nentially decreasing hopping is essential to introduce all 
the energy scales while still keeping the numerical effort 
reasonable. Such a procedure is especially suitable for 
the Kondo problem where a range of energy scales con- 
tributes equally to the screening of the impurity spin. 
The algorithm is stopped after iVmax iterations. In the 
presented results we have typically used A = 3 — 4, 
K = 2000 (not counting the degeneracies due to spin, 
isospin, and parity symmetries'^*' which have been explic- 
itly taken into account) and iVmax = 40. 

In order to gain additional insight and to make a rela- 
tion with our previous work we compare the NRG results 
to the results obtained by the Schonhammer-Gunnarsson 
(SG)^i^ variational method. The details of our imple- 
mentation of the variational method are given in our 
previous work^Si^iM. For reader's convenience we here 
just remark that it consists of finding the parameters of 
an auxiliary noninteracting Hamiltonian H [of the same 
form as H in Eq. (jT]), but ior g — 0,U = and described 
by renormalized parameters VL,VR,e\, which minimize 
the variational ground state energy E = (^|if |\1/), where 
the variational function ^ is expressed in the basis of pro- 
jection operators Pi acting on the Hartrec-Fock ground 
state 1^0) (which includes the phonon vacuum) of the 
auxiliary Hamiltonian H, 



\^)=Y.^n 



(an"P.|*c 



(14) 



We have adapted the SG method also to extract the 
effective oscillator potential. By restricting the param- 



For U and fi large the low energy behaviour is obtained 
by projecting the Hamiltonian onto space consisting of 
states with singly occupied molecular orbital and without 
excited phonons. The result of this Schricffer- Wolff (SW) 
transformation (described in Appendix B) is the 2CK 
Hamiltonian 



H- 



2CK 



— Je3 • Sg -\- JqO • So, 



(15) 



describing the anti-ferromagnetic coupling, between the 
spin on the molecular orbital S and the spin densities s^ 
in orbitals next to impurity in the even and odd channel, 
a = e, o, respectively. The coupling constants are 



J, = 2K^ ( — + — 
-e e + 



U 



and 



Jo = 2K' 



1 



1 



-e + Q e + U + n 



The ratio between the coupling constants 
Jo ^ g^ 

Je 1 -I- 2niu 



(16) 



(17) 



(18) 



determines in which of the two channels the Kondo 
screening takes place. Let g^ denote the delimiting 
value separating regimes with different symmetries of 
the screening channel. At g^ both channels participate 
equally to the screening and the over-screened non- Fermi 
liquid behavior results. In terms of the original model, 
this corresponds to the point at which inequality Eq. ^ 
is violated. According to Eq. (fT^ . for large C/, J7, 



9c 



yji + 2n/u. 



(19) 



We now turn to the renormalization of the vibrational 
mode and demonstrate that through the electron-phonon 
coupling the confining potential is diminished and can 
even be driven to the form of a double well. We first 
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Figure 3: (Color online) (a) Effective oscillator potential for 
LM. Semi-classical estimate (dotted); SG estimate for C = 
(full) and C = 0.01 (dashed). Parameters f7 = 0.1, V = 0.02 
and U — 0.3 (for SG only) are in units of D (half-width of 
the band), (b) The minima of the potential in semi-classical 
estimate. 



discuss U ~ model and use the semi-classical ap- 
proximation in which phonon operators are substituted, 
a ~ flt — > x/2, by a real valued constant x. The resulting 
Hamiltonian with the hybridization r(a;) = r(l + g^x'^) 
is noninteracting and can be thus solved exactly. Here 
the magnitude of the bare hybridization with 5 = is 
defined by F = nV'^/D. In this model only the hy- 
bridization energy gain, which in the wide-band limit 
{T/D smah) reads?-?^ AE^yh = -2/7rflogL»/f, and the 
elastic energy cost AE^i — Vtx"^ /A are a function of x. 
Hence, the effective oscillator potential in this estimate 
is AEgc = AEc\ + Ai^hyb and can be written in a closed 
form 

A£;sc(a;) = nx^A - i2/n)t{x) log{D/[fix)]}. (20) 

The prefactor of the x^-tcrm in the snrall-a; expansion is 
equal to Q/A- {{2/n)g^T [log(D/F) - 1]} and is decreas- 
ing with increasing g indicating the softening of the con- 
fining potential. At sufficiently large g, two wells emerge. 
In Fig. [3] (a) we plot (dotted) the resulting potential for 
g below and above the delimiting value (gd)- The gd and 
the positions Xm [shown in Fig. [3{b)] of the potential 
minima can be extracted analytically from Eq. (j20p . 
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9d 



8F [logiD/T) - 1] ■ 



Un{g- gd) 



4r<?^ 



(21) 
(22) 



For a finite U the model cannot be solved exactly, but 
by estimating the potential with the SG method, we find 
that the evolution of the potential as just described per- 
sists. The SG results for [/ > are plotted in Fig. [3] for 
C, — (full lines) and also for the case with broken left- 
right symmetry C, = 0.01 (dashed). While for g < gd the 
potential is only slightly perturbed, finite C for g > gd 
breaks the degeneracy between the two minima. In this 
regime, the molecule is attracted to one of the electrodes. 



Table I: SC and SG: (?d obtained from semi-classical (SC) es- 
timate and from SG simulations. SW and NRG: Qc calculated 
by SW estimates and from NRG simulations. SG method is 
inapplicable for very large U. 
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2.76 
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2.60 
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1.29 


1.28 



Having defined the characteristic values gc and gd it is 
interesting to ask whether there is any relation between 
the two. As discussed above, the 2CK point occurs at gc 
such that in the semi-classical description 1 — gda^l ~ 0, 
meaning that the double well-potential has to be pre- 
formed for its minima to occur at Ixml ^ l/5c- Therefore 
we expect gc > gd- However, as Xm evolves rapidly as a 
function of g the values of gc and gd are close. 

To support this statements quantitatively in Table H] 
we show the semi-classical (SC) and SG estimates of gd 
compared to the SW and NRG estimates of gc for sev- 
eral values of parameters. As expected, we find that 
the Schrieffer- Wolff estimate of gc becomes more accu- 
rate (agrees better with the NRG) for U, fl large when 
the charge fluctuations and phonon excitations are sup- 
pressed. Conversely, the semi-classical estimate is more 
accurate for large number of excited phonons (small fi) 
and small U but breaks down for large U. For exam- 
ple, for fi = 1 and U = 3 the value obtained by the 
semi-classical method overestimates gd to a value which 
is larger than gc obtained using NRG. 



B. Exponential modulation 

Let us perform the Schrieffer- Wolff transformation also 
on the model with exponential modulation of tunneling. 
We obtain 



Je,o = "iV 



00 

m=0 



-e + raVl e + t/ + raVl 



(23) 



where (5^ = |(0| cosh(5a;)|m)p and (5°„ = 
|(0| sinh(5a;)|m)p. Note that as Je > Jo for all 
values of parameters no 2CK fixed point occurs in 
such a model. Note also that both J^ depend on g 
exponentially. The Kondo temperature, which itself 
is exponential in J, Tk oc exp(— 1/pJ), is thus very 
sensitive to the value of g. 

The effective potential which is unbounded from be- 
low for L — > cxD is regularized by the phonon cutoff L. 
At fixed L we define gdL) as the value of g at which 
the molecule is attracted by the hard wall boundary as 
signaled by an abrupt increase of displacement (shown 



later). Given L, the behavior of the model for g > gdL) 
becomes doniinated by the hard-wall boundary only. For 
instance, the full curve in Fig. [2Ub) corresponds to a par- 
ticular g which is smaller than g^ for the smaller phonon 
cutoff (dotted region) and larger than g^ for the larger 
phonon cutoff (dashed region). For i — > oo, gc ^ 0. 



V. NUMERICAL RESULTS 

In this section we confirm the anticipations stated 
above with numerical examples. We first show the re- 
sults for LM and then for EM. We treat separately the 
inversion symmetric C = ^^^d asymmetric C > cases. 
We use the half-width of the band D as the energy 
unit. Unless where explicitly stated, we take U = 0.3, 
r = nV^/D = 0.02, and Q = 0.1. All the results corre- 
spond to the particle-hole symmetric point, e = — C//2, 
and to the zero temperature limit. 
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Figure 4: (Color online) (a) Fluctuations of displacement, 
(b) Fluctuations of charge, hopping to one of the electrodes 
and its fluctuations. The left (right) pictograms schemati- 
cally present the efi'ective oscillator potential before (after) 
the emergence of the soft mode {U — 0.3, F — 0.02, fl = 0.1). 



A. Linearized model 

1. Inversion symmetry: ^ ~ 

We begin by looking at the static quantities for C = 0. 
The average displacement (x) which is odd under inver- 
sion vanishes. The fluctuations of displacement Ax — 
{{x — (x))^)^/^, shown in Fig. [H^a), increase monotoni- 
cally with g. The slope of Ax increases considerably at 
g ^ g^ (~ gj.)^ where the double well effective poten- 
tial is formed, as indicated in pictograms. The change 
in slope is driven by the increased hybridization in the 
odd-channel. The position of the minima of effective po- 
tential Xm (full line) Eq. (|22p also becomes non-vanishing 
there. 

Due to increased hybridization the fluctuations of 
charge An^ = ((n— 1)^), shown in Fig.[3]^b), are increased 
in the g > gc regime. However, the absolute value of the 
average of the hybridization operator vl = (vl) = (vr) 
is diminished, contrary to the naive expectation. This is 
due to increasingly fluctuating sign of the overlap integral 
in this regime. On the other hand, the average of the hy- 
bridization operator squared w£ — (^\) is increased here, 
as expected. 

In Fig. [5] we plot the NRG flow diagram: the energies 
of the lowest few eigenstates in units of characteristic 
energy of a particular iteration ujj^ (x K'~^'^~^^/'^ as the 
function of the NRG iteration number N . The finger- 
print of the Fermi liquid ground state are the equidistant 
low-lying quasiparticle excitation^, which are seen for 
large N , irrespective of g. By comparing the top two 
panels with the bottom panel it is seen, that the roles 
of even and odd parity states are interchanged in the 
Fermi-liquid regime (right-hand side of each panel) cor- 
responding to the change of the screening channel as g is 
increased above g^- 

For g ^ gc (bottom panels) the unstable non-Fermi liq- 
uid fixed point, which determines the NRG flow at inter- 
mediate N (~ 20 for the plotted case) is discerned. The 
ratios of eigen-energies here (horizontal bars) are charac- 
teristic of the 2CK effect and agree with the predictions of 
conformal field theory^e (o, 1/8, 1/2, 5/8, 1, 9/8, ...). This 
regime cannot be explained in terms of the Fermi liquid 
quasiparticles. The difference between the couplings to 
the screening channels is a relevant perturbation and for 
low temperatures (large TV) drives the flow towards the 
(stable) Fermi liquid fixed point. 

Now we turn to the renormalization of the phonon 
propagator by the electron-phonon coupling. The dy- 
namical information about oscillator is contained in the 
displacement Green's function. The displacement spec- 
tral function 

A[iLJ) = Im <^ x,x >c^ = 



--Im 



{^{){[x{t),x{Q)])e'^'dt (24) 



is an odd function of lu due to the hermiticity of x (unlike 
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Figure 5: (Color online) NRG flow diagram. In the left (right) 
panels eigen-energies of states of even (odd) parity are shown 
for odd (even) number of NRG iterations, respectively. The 
states are labeled by {Q,S): total charge Q and total spin S 
quantum numbers. The thick horizontal bars are the confor- 
mal field theory predictions for 2CK fixed point. Parameters 
are as in Fig. |4] 



Ini '^ a,a'' 3>[j which is odd only for the inversion sym- 
metric case C = 0). Since in NRG A{uj) is evaluated for 
a finite system it consists of several (5-peaks of different 
weights. To obtain a smooth spectral function we have 
used the Gaussian broadening on the logarithmic scaled, 
where the Dirac S function is broadened according to 
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, (25) 



and we used b — 0.3 in our calculations. 

In Fig. [HJa) we plot A{oj) for various g. The width 
of the high frequency peaks is overestimated (the ex- 
treme example is the g = peak at f2 for which the 
width should vanish) due to the broadening procedure 
described above. We could use Dyson equation^i^S to 
obtain sharper peaks but on one hand there is no a ■pri- 
ori guarantee that such a procedure gives more accurate 
results for large g and on the other hand we do not use 
the width of the peaks as a means to draw any quantita- 
tive conclusion. 

The evolution of the phonon operator can be under- 
stood in terms of the evolution of effective confining po- 
tential (Fig.[3l pictograms in Fig. |4]). For intermediate g 
[starting at g ~ 0.5 for the parameters used in Fig. EJa)] 
as the confining potential starts to diminish the vibra- 
tional mode begins to soften: the peak of A{Ld) moves 
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Figure 6: (Color online) (a) Displacement spectral functions 
for various g. Note the softening of the phonon mode. For 
large g only small amount of spectral weight resides at low 
frequencies. The peaks at large frequencies appear broader 
than they should be because of the broadening procedure de- 
scribed in the text. Parameters are as in Fig. |4l (b) Spectral 
function after the emergence of the soft mode schematically. 
The frequencies of the soft mode uiq and the high-frequency 
oscillation at oj ~ SI are indicated. 



to lower frequencies. At still larger g ^ gd two peaks 
emerge at the point where the two wells develop in the 
effective potential. Characteristic dependence of A{uj) in 
this regime is schematically presented in Fig. El^b). The 
double well potential is already well established and the 
major part of the spectral weight corresponds to the os- 
cillations within each of the wells. The minor, low fre- 
quency part corresponds to increasingly slow tunneling 
(see also Fig. [8]) between the degenerate minima of the 
potential. The frequency wq and the weight of the low 
frequency peak decrease with increasing g. 

Due to the inversion symmetry the conductance can 
be calculated using the scattering phase shifts only, as 
described in Appendix C. The scattering phase shifts in 
the even and odd channels are 7r/2(0) and 0(7r/2), for g < 
5c {g > 5c) and the conductance evaluated by Eq. (jC8|) 
is unity2i. 



2. Broken inversion symmetry: ( > 

It is impossible to experimentally produce perfectly 
symmetric devices, therefore it is interesting to check for 
the influence of the inversion symmetry breaking term of 
relative strength (. Let us first remark that the NRG 
fiow diagrams (not shown here) arc that of the Fermi liq- 
uid as the occurrence of 2CK fixed point is inhibited by 
the breaking of inversion symmetry. 

In Fig. [7] we plot static correlations for ( — 0.01. 
New compared to the inversion symmetric case is the 
non- vanishing average displacement, which monotoni- 
cally increases with g. Despite the simultaneous in- 
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Figure 7: (Color online) Static quantities for system with bro- 
ken inversion symmetry, (a) Averages of x, x^ , and Aa:. (b) 
Averages of hopping vl and vr and fluctuations of charge. 
The breaking of inversion symmetry C, = 0.01. Other param- 
eters as in Fig. |4] 



crease in average of x^ , the fluctuations of displacement 
Aa; = {{x — (x))^)^'^ eventually reach a maximum. At a 
still larger g the molecule is attracted to the right elec- 
trode as indicated in the right pictogram. The fluctu- 
ations of charge are remain the same as in the ^ = 
case but the average hopping to the right electrode vji 
is larger than v^, which first vanishes and then changes 
sign. This asymmetry happens not at gc ~ gd but at 
another value gasy {'^ 1 for these parameters), when the 
softened phonon frequency ujq and the energy difference 
between the hybridizations to left and right ex C become 
comparable. 

The spectral functions for ^ — 0.01 differ from the 
C = case only for g > gasy The distinction between 
the two cases is mainly that for C = 0.01 the frequency 
of the soft mode oscillation ojo saturates. In Fig. [SI we 
plot the ujQ for both ^ as a function of g. For ( — 0.01, 
we also plot the weight of the luq peak, which diminishes 
exponentially with increasing g. 

Again, the behavior is easily understood in terms of the 
effective potential pictorially shown in Fig. \7\ For small 
to intermediate g the soft mode begins to emerge as the 
shape of the potential is transformed to a double-well- 
like form with the right well being lower in energy by a 
value A (x C. When g increases further and loq decreases 




Figure 8: (Color online) The frequency (full lines) of the soft 
mode peak as a function of g. The weight of the soft mode 
peak for (^ = 0.01 and normalized to some arbitrary value 
(dashed) . Other parameters are as in Fig. 3] 







Figure 9: (Color online) (a) Displacement spectral functions 
for flxed g = 1.05 and ( = 0,0.02,0.04. Note the reduction 
of weight of low frequency peak on increasing C,. The high 
energy peak of ^ = case is approached at still smaller C,. 
Inset: The position of low frequency peak as a function of 
C for n = 0.1, g = 1.05 and STJ = 1, 5 = 7. The fuU lines 
are fits to the data with slopes 0.32 and 0.42 for fi = 0.1, 1, 
respectively. Other parameters are as in Fig. 21 



below A, the tunneling is suppressed. In this regime, the 
major part of the displacement spectral weight is due to 
the tunneling within the lower of the wells. The average 
displacement increases, its fluctuations decrease. Note 
that even a minor breaking of inversion symmetry results 
in a strongly asymmetric state through the mechanism 
described here. 

We plot also A{lo) for different C, and fixed g = 1.05 > 
5asy in Fig. HI The characteristic frequency of the soft 
mode is related to the energy-difference of the two wells 
and is proportional to C, as shown in the inset of Fig. [9l 

The conductance for C, > Q cannot be obtained from 
the scattering phase shift alone since the rotation angle 6 
(see Appendix C) is not known in the present case of bro- 
ken symmetry. Therefore we calculated the conductance 
from the current-current correlation function. We plot 
the conductance calculated by NRG and SG in Fig. [TOl 
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Figure 10: (Color online) Conductance as calculated from the 
current-current correlation function obtained by NRG (full) 
and from the SG method. The shaded region indicates the 
region the results are outside of the scope of the linearized 
model. C, = 0.01, other parameters are as in Fig. 3] 
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The conductance decreases to zero at g ^ gc- This min- 
imum corresponds to the virtual decouphng of the left 
electrode, 1 — gx ^ Q, there on the average. Due to lin- 
earization, for still larger g the magnitude of the overlap 
to the left electrode increases again. Correspondingly, 
the conductance is increased. The access to this region 
(denoted shaded) is probably unaccomplishable in mea- 
surements of transport through molecules. 

The NRG and the SG data agree well for most g. How- 
ever, for g ~ (7asy the maximum and minimum are seen in 
the NRG results, a feature which the effective Hamilto- 
nian of the SG method fails to capture. The discrepancy 
is especially visible for the parameters used here for which 
5c ~ 5asy (compare also with data given in Appendix A) . 
The NRG and SG for 5 > ^asy agree again. Here the 
fluctuations of displacement are diminished and the be- 
havior is efficiently described in terms of the effective 
Hamiltonian with asymmetric coupling to the electrodes. 



B. Exponential model 

We now turn to the model with exponential modu- 
lation. Here even for the inversion symmetric model 
(C = 0) the 2CK Kondo fixed point is inaccessible because 
the coupling to the even channel is invariably stronger 
than the coupling to the odd channel. Correspondingly, 
the finite size spectra shown in Fig. [TlTa) are that of the 
Fermi liquid ground state with the spin screened by the 
even conduction channel. Even at large g, the even and 
the odd channels do not inter-change their roles in the 
screening, confirming such an intcr-change in the model 
with LM is indeed an artifact of the linearization. 

In Figs.fTlTb) we plot the impurity contribution to the 
magnetic susceptibility x- Dimensionless susceptibility 
A:Tx/(2^b)^, where k is the Boltzmann constant and jjlb 
the Bohr magneton, has a peak at intermediate temper- 
atures corresponding to the local moment regime pro- 
vided that g is below some critical value g^- For g > g^ 



Figure 11: (Color online) Model with exponential modulation 
of overlap integrals (EM). L — 10, other parameters as in 
Fig. [l] (a) NRG eigenvalues, (b) Impurity contribution to 
magnetic susceptibility. For g > gc the local moment regime is 
absent. Diamonds indicate the susceptibility calculated using 
the Bethe ansatz. (c) The Kondo temperature. 



the local moment regime is absent. By fitting the Bethe 
ansatz results for S = 1/2 Kondo impurity to the numer- 
ically calculated susceptibility, following the procedure 
described in Refs. l4Clll4lll43 . we obtain the estimate of the 
Kondo temperature Tk, shown in Fig. Illf c). Tr- is in- 
creasing rapidly near gc as effective hybridization grows 
large. For g > g^ the Kondo temperature is not defined 
because there is no local moment in the system. 

In Fig. [iWa) we plot the averages of displacement and 
its fluctuations for C = 0.1. For g > gc, the displacement 
rises abruptly to the values that arc limited only by the 
phonon cutoff; for C — ^ (i^ot shown) the same applies 
to the fluctuations of displacement (displacement itself is 
zero). The abrupt increase is due to the increased hy- 
bridization, which is exponential in displacement. The 
gain in the kinetic energy cannot be compensated by 
the cost in oscillator potential, which is only quadratic 
in the displacement operators. In this regime, due to 
the exponentially increased hybridization, the fluctua- 
tions of charge, shown in Fig. I12r b). are near-maximal. 
The molecule resides in the effective potential presented 
in the right pictogram, which has the same form as the 
function plotted in Fig. [D^b) (dashed), and is attracted 
to the right electrode. 

The value gc can be estimated by first recognizing that 
the maximum x is bounded by the phonon cutoff L: a; < 
a;max = '^VL and then solving flxf^^jA = V^exp(ga::max) 
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Figure 12: (Color online) (a) Average displacement and its 
fluctuations, (b) Fluctuations of charge. Phonon cutoff L — 
10, (" = 0.1, other parameters as in Fig. |4l 



for g, which gives 



^° = ^^°K^'^ 



(26) 



The comparison between this estimate and the numerical 
data is shown in Fig. [T31 In the limit of large number of 
allowed phonons L -^ oo critical value g^ —^ 0. An im- 
portant result is that the model with exponential modu- 
lation of hybridization is not well-defined with quadratic 
stabilizing potential only as the results strongly depend 
on the cutoff. 

In Fig. [T3]we plot the spectral function A{uj). With 
increasing g ^ gc the bare oscillator peak starts to soften. 
At g > gc the molecule is attracted next to the hard wall 
boundary and remains mainly trapped into one of the 
wells defined by the phonon-cutoff. The high-frequency 
part of the spectral function in Fig. [14] corresponds to 
the oscillations within these wells. The wells are strongly 
anharmonic which is reflected in the broad distribution 
of spectral weight. The low frequency part of the spectral 
weight is due to the oscillations between the wells. 

A natural question which arises at this point is whether 
it is possible to tune the parameters so as to drive EM to 
the regime with developed double well potential, but for 
g sufficiently lower than g^, so that the wells are not next 
to the hard wall boundary. For the parameter set used 
in the results shown, for example, this is not possible, 
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Figure 13: (Color online) Values ^c where the breakdown 
of the exponential model occurs for phonon cutoffs L = 
10, 20, 30, 40, 70 (circles). Other parameters are as in Fig. |4l 
Semi-classical estimate of g^ (dashed). 
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Figure 14: (Color online) Spectral functions exponential mod- 
ulation. (" = 0.1, L = 10, other parameters as in Fig. [l] 



because 



5d > 5o 



(27) 



for all L. 

It can be shown that the inequality Eq. ((27)) holds in 
general. We begin by maximizing g^ with respect to L 
(treating L as if it was a continuous variable) and obtain 



g, = max.9e(L) = ^^^j;^. 



(28) 



The ratio between g^ defined in Eq. (|2ip and g'^ can now 
be evaluated. We obtain that gd/g'^ > (e'^7r)i/'*/4 - 1.92. 
The double well potential thus develops also in EM but 
the wells are near the boundary defined by the phonon 
cutoff. 



C. Regularized exponential model 

By using EM some of the nonphysical results found in 
the model with LM are eliminated but others, such as 
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Figure 15: (Color online) Results for a model with over- 
lap integrals modulated as Vl,r{x) ~ exp(=p(;3;)/cosli((7a;), 
Eq. (|29p . with ( — 0.1 and other parameters as in Fig. |4] (a) 
Displacement and its fluctuations, (b) Expectation values of 
hoppings and fluctuations of charge. 



the dependence on a cutoff parameter L, are introduced. 
Another possibility is to regularize EM, 

VlMx) = V [expiTgx)/ cosh((7x) T C] , (29) 
or in the symmetrized basis 

Ve ^ V2V, Vo = V2V [tanh(5a;) + C] • (30) 

The inequahty Eq. ^ is satisfied and the normahzation 
with the cosh function ensures the model behaves well 
for large x. 

For this model we evaluated the matrix elements of the 
Hamiltonian {m\H\n) {\'m), \n) correspond to states with 
TO, n excited phonons) in the real space by reintroducing 
the Hermite polynomials, which is numerically more sta- 
ble than the expansion of tanh(gx) in the power series 
in X. This procedure can be used for any form of modu- 
lation. For EM the procedure can be simplified because 
it is possible to evaluate (m| exp((7a;)|n) analytically via 
the Baker-Hausdorff equality. 

We first present results for parameters kept as in LM, 
Fig. |4l and with asymmetry parameter <^ = 0.1. Due to 
weaker dependence of the overlap integrals on displace- 
ment the effects of the electron-phonon coupling in this 
model are less pronounced. The displacement and its 
fluctuations, which we plot in Fig. [TsT a) are small. This 
is accompanied by only a minor softening of the phonon 
mode, as shown in Fig. [TCT a). The fluctuations of the 
charge and the expectation values of hopping, plotted in 



Figure 16: (Color online) (a) Displacement spectral functions 
for standard parameters as in Fig. [K] (b) Spectral functions 
for the same set of parameters, except for softer spring con- 
stant Q. = 0.01 and C. = 0.01. 



Fig. [TST b) are likewise only minorly distorted from the 
g — Q case. Note that we used C = 0.1 here, therefore 
the hoppings toward left and right electrodes differ con- 
siderably already for the g = Q case, a feature which is 
only slightly (compared to LM and EM) amplified by the 
dependence of the overlap integrals on x for g > 0. The 
reason for this moderate dependence of quantities on g is 
the fact that the largest energy the system can gain by 
increasing the displacement is of the order of V oc vT, 
which is for these parameters comparable to the elastic 
energy Oa;^ already for the displacements of order 1. 

For smaller 17 = 0.01, i.e., a softer 'spring', the effect of 
the electron-phonon coupling is larger and the soft mode 
is clearly developed. Fig. [TBT b). Simultaneously other 
quantities ((a;), Ax, etc.) also exhibit more pronounced 
behavior, similar to LM and EM cases (not shown here). 

The regularized form, Eq. [291 might describe well the 
modulation of overlap integrals in realistic case and we 
believe that the softening of the phonon mode and re- 
lated tendency toward broken symmetry configuration 
is the universal consequence of displacement-modulated 
hybridization. 



VI. CONCLUSION 

In summary, we analyzed the influence of the electron- 
phonon coupling in molecular bridges consisting of a 
molecule oscillating between two electrodes. The over- 
laps between the molecular orbital and the orbitals in the 
electrodes are determined by the position of the molecule 
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X. To model this situation we used several types of the 
dependence of the overlap integrals on x. 

We find that the inversion symmetric model with lin- 
ear modulation has a 2CK critical point at some criti- 
cal electron-phonon coupling g — gc where both chan- 
nels participate equally to the screening of the spin. The 
occurrence of this critical point is suppressed if a finite 
difference C between the coupling to left and right elec- 
trodes is introduced. In such broken symmetry system 
the sharp transition between the two Fermi liquid states 
via a non Fermi liquid state is replaced by a continuous 
rotation of the screening state in the channel space from 
an almost symmetric to an almost antisymmetric linear 
combination of the left and right states. This continuous 
rotation is accompanied by a dip in conductance near the 
point where one of the electrodes is decoupled. 

Additionally we found that the electron-phonon cou- 
pling modifies the shape of the effective potential affect- 
ing the static and dynamic properties of the oscillator. At 
moderate g the potential is softened and the frequency of 
oscillations decreases. At large g the potential develops 
side- wells and the phonon propagator consists of a part 
corresponding to high frequency oscillations within the 
well and another part corresponding to slow tunneling 
between the degenerate wells. 

For finite C, the degeneracy between the wells is broken 
by an amount ex C and when the softened frequency wo 
drops below this value, the tunneling to the higher well 
is suppressed. In this regime, the average displacement 
starts to increase significantly and its fluctuations de- 
crease. Hence, only a minor (~ wq) breaking of inversion 
symmetry can result in a significantly asymmetric con- 
figuration, which could also account for the asymmetric 
configurations typically observed in some experiments. 

We consider also overlap integrals exponentially modu- 
lated by the displacement, which is closer to reality in the 
respect that the coupling to the even channel is invariably 
stronger. However, due to the exponential increase in hy- 
bridization energy gain, this model is not stable against 
large distortions. The value of g at which this instabil- 
ity starts is set by the phonon cutoff and vanishes in the 
limit of large cutoff. We analyzed the model at finite 
cutoff and found that the 2CK fixed point is absent, but 
the other behavior of the LM is qualitatively recovered. 
In particular, the softening occurs with increasing g and 
for large g the molecule is attracted to one of the elec- 
trodes and is localized next to the boundary given by the 
phonon cutoff. In this regime, the Kondo temperature is 
significantly increased but the conductance is suppressed 
due to the small overlap with one of the electrodes. 

The main common finding is thus the softening of the 
phonon mode, emerging for linearized and exponential 
modulations, as well as for the case of a more realis- 
tic regularized modulation. If the inversion symmetry 
is broken, the instability due to the formation of a dou- 
ble well effective potential will manifest as the attraction 
of the molecule into one of the wells and simultaneous 
suppression of the conductance. The instability inspired 



also a very recent work, in which the break-junctions are 
studied as an example of a two-level system™. 

Finally, let us comment on the relevance of the pre- 
sented work outside the scope of the experiments with 
molecular conductors. Within the dynamical mean field 
theory (DMFT)^^ the bulk correlated electron systems 
are solved by mapping onto impurity problems. Likewise, 
bulk systems with electron-phonon coupling are mapped 
onto impurity (or impurity-cluster) problems with cou- 
pling to phonons. In this regard the knowledge of the 
behavior of the impurity problems is a convenient guide 
in the interpretation of the DMFT results. The results 
obtained in this work for the linearized model correspond 
to the general two-band case, where the coupling to one 
of the bands is phonon-assisted. The large g regime which 
we dismissed as unphysical could prove relevant in this 
context. 

Our results could also be applied to the studies of 
nanoelectromechanical systems^^'^i^'^i^ (NEMS). In 
NEMS the tunneling to electrodes is modulated by the 
displacement of the cantilever in a similar fashion as ana- 
lyzed in this work. Once the dimensions of these devices 
are reduced to such an extent that the frequencies of the 
oscillations will become comparable to other scales, such 
as the bias at which the devices are operated, softening 
of the vibrational mode and susceptibility towards large 
displacements could be observed. For instance, we pre- 
dict that the frequency of the oscillations will decrease 
if the tunneling rate from the electrodes to the perpen- 
dicularly situated cantilever immersed between them is 
increased. 
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Appendix A: COMPARISON TO 

SCHONHAMMER-GUNNARSSON 

PROJECTION-OPERATOR METHOD 

In this appendix we compare the results of NRG cal- 
culations to the results obtained by the SG method. Let 
us summarize first the results we have obtained using the 
SG method for model II of Refs. [Hill. We have found 
that for g large enough the variational solution with bro- 
ken inversion symmetry is lower in energy. This occurs 
even for C = when the inversion symmetry should per- 
sist. That indicates that for large g due to the instability 
in the system the SG method fails giving a solution with 
'spontaneously' broken symmetry. Such failure is typical 
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Figure 17: (Color online) Comparison between Shonhammer- 
Gunnarsson and NRG results, (a) Displacement and displace- 
ment fluctuations, (b) Conductance and charge fluctuations. 
U = 0.3, r = 0.02, Q = 0.1, C = 0.01. 



of the mean- field treatment. However, from our previous 

analysis^^. it was not clear whether the failure occurs at 
the point where the symmetry of the screening channel 
is changed or at the point when the soft mode is formed 
(or the two phenomena occur simultaneously). 

In Fig. [T7] we compare the displacement, its fluctu- 
ations and the fluctuations of charge calculated by SG 
method to the NRG results. They match closely, with 
the exception of discrepancies in the precise values of in- 
creased displacement fluctuations and range of g where 
these occur. Conductance is discussed in the main text, 
here we re-plot the curve for completeness. 

In Fig. [TST a) we show the results for fi = 1 where the 
discrepancy between NRG results (thick) and SG results 
(thin) is larger. In this regime SG method overestimates 
the value of displacement and its fluctuations. More in- 
terestingly, in Fig. fTSTb) the jump in expectation values 
of hopping operators and a minimum in the conductance 
are seen in the variational results near the value Qc ~ 2.3. 
Here also a small peak in displacement is seen. 

If ^ = there is no 'spontaneous' symmetry breaking 
in SG method for this fi. By combining these results we 
conclude that the appearance of asymmetric solution in 
SG will coincide with the change in the symmetry of the 
screening channel from even to odd, but only when the 
soft mode is sufficiently developed, otherwise only a flnite 
jump in vl, vr occurs there. 

Let us remark here also that in terms of the effective 
Hamiltonian, the g < gc{g > gc) regimes correspond 
to the state, where the hopping to left and right elec- 




Figure 18: (Color online) f7 = 1, other data as in Fig. 1171 
Thick lines: NRG, thin lines: SG; (a) Expectation values of 
X and its fluctuations, (b) Fluctuations of occupancy and 
expectation values of hoppings to left and right electrodes. 
The vertical lines denote gc (thick) and the value of g where 
the minimum of conductance calculated by SG method occurs 
(dotted). Conductance curves (dashed) obtained by NRG and 
SG coincide. 



trodes has equal and opposite phases, respectively. The 
non-Fermi liquid g = gc regime cannot be described in 
terms of (Fermi-liquid) effective Hamiltonian. Neverthe- 
less, insisting upon this description, it can be regarded 
as a combination of states which correspond to two effec- 
tive Hamiltonians in which left and right electrodes are, 
respectively, decoupled. 



Appendix B: SCHRIEFFER- WOLFF 
TRANSFORMATION 



To obtain the effective-low energy Hamiltonian iJeff we 
flrst divide the Hamiltonian into two parts^ 



H = Ha + XH' 
where Hq, which for our example reads 

Hq = en + Un-^ni + ila'a + Hl + Hn, 



(Bl) 



(B2) 



is diagonal in the low {n = 1, no excited phonons) and 
high {n — 2, 0, with excited phonons) energy subspaces 
(for e ^ —U/2, U and fi large). The hybridization part 
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H' = Vve + Vgxvo 



(B3) 



provides the mixing between the low and high energy 
subspaces and A serves as an expansion parameter to be 
set to 1 at the end of the derivation. 

Then a canonical transformation generated by some 
unitary operator e is performed to obtain the block- 
diagonal 



H = e^He^^ = 



Hl 
, Hh. 



(B4) 



By expanding Eq. (|B4p in terms of nested commutators 
H = H+[S,H] + [S, [S, H]]/2 + ... and the generator S 
as a power-series in A, 5 = SiX + S2\^ + ■•■, to second 
order in A the following should hold: 



H = Ha + \{H' + [Si,Ha\)- 



(B5) 



\WSi,H'] + -[Si, [5i,ifo]] + [S2,Ho]). 



Relation Eq. (jB5P is satisfied to lowest order in A by 
demanding 



H' + [SuHq]=Q. 



(B6) 



Since Si can be chosen completely block-off-diagonal the 
commutator [Si, H'] is block-diagonal, therefore Eq. (|B4p 
is satisfied also to second order in A by taking ^2 = 0. 
Looking now at the matrix element of Eq. (|B6p in the 
basis of eigenstates of fl^o (we will use Latin indices 
for states in low-energy subspace, and Greek indices 
for states in high-energy subspace), we obtain Sab = 
H'ab/iEa - Eb). 

Finally, the high energy part is neglected by projecting 
to the low energy subspace 

^off = PiowHqPiow + HeS (^7) 

with the projector 

Plow = |0)(0) K(l - nx) + nx(l - nt)] (B8) 

where [m) denotes normalized phonon state with m ex- 
cited phonons and the effective interaction H\^ is due 
to the virtual transitions to the high energy states. Its 
matrix elements read 



7ehigh \ 



rr/f U' Tin ff' 

— E'y Eb — -C/-7 



(B9) 



In this case the Hamiltonian can be up to a constant 
term recast by using the spin operators to the form of 
the 2CK model 



HcS — 'f^2CK — Je^ ■ Se + Jo^ ' ^o 



(BIO) 



where the Sa for a = e,o denote the spin densities which 
read 



S= -^d\(Tss'ds' 



(Bll) 



for d-orbital and likewise for orbitals a. Here the compo- 
nents of <T are the Pauli matrices. The coupling constant 
to the even channel is 



J, ^2V'\ — 



1 1 



-e e + U 
For the odd channel we get 

2 2^1 1 



Jo = 2V'g 



-e + fi e + t/ + fi 



(B12) 



(B13) 



where the Vi in the denominators occur since high-energy 
states involve one excited phonon. 

For the exponential modulation the hybridization parts 
read 



H'^ =V iiin\i[gx)vo- 



(B14) 



The low-energy Hamiltonian is again that of the 2CK 
model and the coupling constants read 



J. = 21/. 



2 Y^ 



'^-^ V ~e + Tifi ' e + [/ + m^ 



Jo = 2V' 



S° 



SI 



m—O 



, (B15) 



(B16) 



where 6ff = |(0|(e9=^ ± e-9^)/2[m)[^ = exp{g^)g^"^/ml 
for m even (odd) and zero otherwise, respectively. 



Appendix C: CONDUCTANCE 

The linear response conductance for a general interact- 
ing system is given by the Kubo formula!^ 



G= lim —ImCiUuj), 



(CI) 



where 



Cil(w) = -^lim / [/(i),/(0)]e'("+''')*dt (C2) 

'7^0 Jo 

is the retarded current-current correlation function and 
the current operators are defined by the time derivative 
of the electrons in the electrodes 



Nr-Nl 



(C3) 
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For Fermi liquid systems at T = the current-current 
correlation function can be expressed in terms of the 
Green's function G„„' (uj) involving sites n, n' in the elec- 
trodes. The conductance is then given by the Landauer 
formula (Fisher-Lee relation^^) 



G-Go|i(0)p 



(C4) 



where Go = 2e^//i is the quantum of conductance and 
the transmission amplitude t(a;) can be written as^i^ 



tiu) 



1 



-i7Tp{uj) 



"''("-"')G„.„(L^-(-ir;). 



(C5) 



Alternatively, the transmission amplitude can be ex- 
pressed also in terms of the scattering matrix Saa' which 
relates the amplitudes of the outgoing wave in electrode 
a to the amplitude of the incoming wave in the electrode 
a' 



5* = 






(C6) 



where r^ {to) are the respective reflection and transmis- 
sion coefficients. The scattering matrix can be rotated in 
the L-R space to the basis of channels (linear combina- 
tions of left and right states) where it is diagonal^ 



usu-^ = 



ci2S„ 








(C7) 



where U = eKpiiOry) e'xp{iaTz) and r^ are the Pauli ma- 
trices. According to the Landauer formula Eq. (|C4[) 
the zero-temperature conductance is determined by the 
transmission probability i(0) = Srl = ^l, hence we ob- 
tain 



G/Go 



\6a-Sb)sin\20). 



(C8) 



The angle 9 determines the maximal value of conduc- 
tance. In the inversion symmetric case 9 = 7r/4, the 
phase-shifts 5a which occur in the even and odd channels 
can then be extracted from the NRG finite-size spectrai^. 

When the system is not inversion symmetric it is not 
possible to extract 9 since the channel indices are not 
tracked. To evaluate the conductance in this case we 
resort to equation Eq. JCTJ. The current operator / = 
{Nl — Nr)/2 is obtained by commuting A^„ = J2ka "i^kaa 
with the Hamiltonian for a = L, R and reads 



N 



'l(R) = -T^L(fl)(a;)2^(c{^(^)^d- /i.e.), (C9) 



where c|^,o^ creates an electron in the orbital next to 
impurity in the left(right) electrode, respectively. 

Alternatively, the conductance can be calculated also 
from the dependence of the ground state energy on auxil- 
iary magnetic flux^'^'^'^ which is obtained by embed- 
ding the interacting system into auxiliary noninteracting 
ring and threading the ring with magnetic fiux <&. Then 
the conductance can be expressed as2^i^ 



G/Go 



TT E{Tr) - E{0) 
2 A 



(CIO) 



where E(^) is the ground state energy of the auxil- 
iary ring with embedded interacting system and A = 
1/ [p{eF)N] the level spacing in the auxiliary ring with 
the density of states at the Fermi energy p{eF) which 
consists of TV sites. This approach is useful especially 
in connection to variational approaches [such as SG, al- 
though in SG we can extract the conductance also from 
Eq. (IC4p directly and obtain same results], where the 
ground state energy is the most reliable quantity. 
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